home *** CD-ROM | disk | FTP | other *** search
/ CU Amiga Super CD-ROM 21 / CU Amiga Magazine's Super CD-ROM 21 (1998)(EMAP Images)(GB)[!][issue 1998-04].iso / CUCD / Programming / Python-1.4 / Python1.4_Source / Modules / cmathmodule.c < prev    next >
C/C++ Source or Header  |  1996-12-15  |  6KB  |  316 lines
  1. /* Complex math module */
  2.  
  3. /* much code borrowed from mathmodule.c */
  4.  
  5. #include "allobjects.h"
  6. #include "complexobject.h"
  7.  
  8. #include <errno.h>
  9.  
  10. #include "mymath.h"
  11. #include "protos/cmathmodule_protos.h"
  12.  
  13. #ifdef i860
  14. /* Cray APP has bogus definition of HUGE_VAL in <math.h> */
  15. #undef HUGE_VAL
  16. #endif
  17.  
  18. #ifdef HUGE_VAL
  19. #define CHECK(x) if (errno != 0) ; \
  20.     else if (-HUGE_VAL <= (x) && (x) <= HUGE_VAL) ; \
  21.     else errno = ERANGE
  22. #else
  23. #define CHECK(x) /* Don't know how to check */
  24. #endif
  25.  
  26. #ifndef M_PI
  27. #define M_PI (3.141592653589793239)
  28. #endif
  29.  
  30. /* First, the C functions that do the real work */
  31.  
  32. /* constants */
  33. static Py_complex c_1 = {1., 0.};
  34. static Py_complex c_half = {0.5, 0.};
  35. static Py_complex c_i = {0., 1.};
  36. static Py_complex c_i2 = {0., 0.5};
  37. static Py_complex c_mi = {0., -1.};
  38. static Py_complex c_pi2 = {M_PI/2., 0.};
  39.  
  40. /* forward declarations */
  41. staticforward Py_complex c_log();
  42. staticforward Py_complex c_prodi();
  43. staticforward Py_complex c_sqrt();
  44.  
  45.  
  46. static Py_complex c_acos(x)
  47.     Py_complex x;
  48. {
  49.     return c_neg(c_prodi(c_log(c_sum(x,c_prod(c_i,
  50.             c_sqrt(c_diff(c_1,c_prod(x,x))))))));
  51. }
  52.  
  53. static Py_complex c_acosh(x)
  54.     Py_complex x;
  55. {
  56.     return c_log(c_sum(x,c_prod(c_i,
  57.             c_sqrt(c_diff(c_1,c_prod(x,x))))));
  58. }
  59.  
  60. static Py_complex c_asin(x)
  61.     Py_complex x;
  62. {
  63.     return c_neg(c_prodi(c_log(c_sum(c_prod(c_i,x),
  64.             c_sqrt(c_diff(c_1,c_prod(x,x)))))));
  65. }
  66.  
  67. static Py_complex c_asinh(x)
  68.     Py_complex x;
  69. {
  70.     return c_neg(c_log(c_diff(c_sqrt(c_sum(c_1,c_prod(x,x))),x)));
  71. }
  72.  
  73. static Py_complex c_atan(x)
  74.     Py_complex x;
  75. {
  76.     return c_prod(c_i2,c_log(c_quot(c_sum(c_i,x),c_diff(c_i,x))));
  77. }
  78.  
  79. static Py_complex c_atanh(x)
  80.     Py_complex x;
  81. {
  82.     return c_prod(c_half,c_log(c_quot(c_sum(c_1,x),c_diff(c_1,x))));
  83. }
  84.  
  85. static Py_complex c_cos(x)
  86.     Py_complex x;
  87. {
  88.     Py_complex r;
  89.     r.real = cos(x.real)*cosh(x.imag);
  90.     r.imag = -sin(x.real)*sinh(x.imag);
  91.     return r;
  92. }
  93.  
  94. static Py_complex c_cosh(x)
  95.     Py_complex x;
  96. {
  97.     Py_complex r;
  98.     r.real = cos(x.imag)*cosh(x.real);
  99.     r.imag = sin(x.imag)*sinh(x.real);
  100.     return r;
  101. }
  102.  
  103. static Py_complex c_exp(x)
  104.     Py_complex x;
  105. {
  106.     Py_complex r;
  107.     double l = exp(x.real);
  108.     r.real = l*cos(x.imag);
  109.     r.imag = l*sin(x.imag);
  110.     return r;
  111. }
  112.  
  113. static Py_complex c_log(x)
  114.     Py_complex x;
  115. {
  116.     Py_complex r;
  117.     double l = hypot(x.real,x.imag);
  118.     r.imag = atan2(x.imag, x.real);
  119.     r.real = log(l);
  120.     return r;
  121. }
  122.  
  123. static Py_complex c_log10(x)
  124.     Py_complex x;
  125. {
  126.     Py_complex r;
  127.     double l = hypot(x.real,x.imag);
  128.     r.imag = atan2(x.imag, x.real)/log(10.);
  129.     r.real = log10(l);
  130.     return r;
  131. }
  132.  
  133. static Py_complex c_prodi(x)
  134.      Py_complex x;
  135. {
  136.     Py_complex r;
  137.     r.real = -x.imag;
  138.     r.imag = x.real;
  139.     return r;
  140. }
  141.  
  142. static Py_complex c_sin(x)
  143.     Py_complex x;
  144. {
  145.     Py_complex r;
  146.     r.real = sin(x.real)*cosh(x.imag);
  147.     r.imag = cos(x.real)*sinh(x.imag);
  148.     return r;
  149. }
  150.  
  151. static Py_complex c_sinh(x)
  152.     Py_complex x;
  153. {
  154.     Py_complex r;
  155.     r.real = cos(x.imag)*sinh(x.real);
  156.     r.imag = sin(x.imag)*cosh(x.real);
  157.     return r;
  158. }
  159.  
  160. static Py_complex c_sqrt(x)
  161.     Py_complex x;
  162. {
  163.     Py_complex r;
  164.     double s,d;
  165.     if (x.real == 0. && x.imag == 0.)
  166.         r = x;
  167.     else {
  168.         s = sqrt(0.5*(fabs(x.real) + hypot(x.real,x.imag)));
  169.         d = 0.5*x.imag/s;
  170.         if (x.real > 0.) {
  171.             r.real = s;
  172.             r.imag = d;
  173.         }
  174.         else if (x.imag >= 0.) {
  175.             r.real = d;
  176.             r.imag = s;
  177.         }
  178.         else {
  179.             r.real = -d;
  180.             r.imag = -s;
  181.         }
  182.     }
  183.     return r;
  184. }
  185.  
  186. static Py_complex c_tan(x)
  187.     Py_complex x;
  188. {
  189.     Py_complex r;
  190.     double sr,cr,shi,chi;
  191.     double rs,is,rc,ic;
  192.     double d;
  193.     sr = sin(x.real);
  194.     cr = cos(x.real);
  195.     shi = sinh(x.imag);
  196.     chi = cosh(x.imag);
  197.     rs = sr*chi;
  198.     is = cr*shi;
  199.     rc = cr*chi;
  200.     ic = -sr*shi;
  201.     d = rc*rc + ic*ic;
  202.     r.real = (rs*rc+is*ic)/d;
  203.     r.imag = (is*rc-rs*ic)/d;
  204.     return r;
  205. }
  206.  
  207. static Py_complex c_tanh(x)
  208.     Py_complex x;
  209. {
  210.     Py_complex r;
  211.     double si,ci,shr,chr;
  212.     double rs,is,rc,ic;
  213.     double d;
  214.     si = sin(x.imag);
  215.     ci = cos(x.imag);
  216.     shr = sinh(x.real);
  217.     chr = cosh(x.real);
  218.     rs = ci*shr;
  219.     is = si*chr;
  220.     rc = ci*chr;
  221.     ic = si*shr;
  222.     d = rc*rc + ic*ic;
  223.     r.real = (rs*rc+is*ic)/d;
  224.     r.imag = (is*rc-rs*ic)/d;
  225.     return r;
  226. }
  227.  
  228.  
  229. /* And now the glue to make them available from Python: */
  230.  
  231. static object *
  232. math_error()
  233. {
  234.     if (errno == EDOM)
  235.         err_setstr(ValueError, "math domain error");
  236.     else if (errno == ERANGE)
  237.         err_setstr(OverflowError, "math range error");
  238.     else
  239.         err_errno(ValueError); /* Unexpected math error */
  240.     return NULL;
  241. }
  242.  
  243. static object *
  244. math_1(args, func)
  245.     object *args;
  246.     Py_complex (*func) FPROTO((Py_complex));
  247. {
  248.     Py_complex x;
  249.     if (!PyArg_ParseTuple(args, "D", &x))
  250.         return NULL;
  251.     errno = 0;
  252.     x = (*func)(x);
  253.     CHECK(x.real);
  254.     CHECK(x.imag);
  255.     if (errno != 0)
  256.         return math_error();
  257.     else
  258.         return newcomplexobject(x);
  259. }
  260.  
  261. #define FUNC1(stubname, func) \
  262.     static object * stubname(self, args) object *self, *args; { \
  263.         return math_1(args, func); \
  264.     }
  265.  
  266. FUNC1(cmath_acos, c_acos)
  267. FUNC1(cmath_acosh, c_acosh)
  268. FUNC1(cmath_asin, c_asin)
  269. FUNC1(cmath_asinh, c_asinh)
  270. FUNC1(cmath_atan, c_atan)
  271. FUNC1(cmath_atanh, c_atanh)
  272. FUNC1(cmath_cos, c_cos)
  273. FUNC1(cmath_cosh, c_cosh)
  274. FUNC1(cmath_exp, c_exp)
  275. FUNC1(cmath_log, c_log)
  276. FUNC1(cmath_log10, c_log10)
  277. FUNC1(cmath_sin, c_sin)
  278. FUNC1(cmath_sinh, c_sinh)
  279. FUNC1(cmath_sqrt, c_sqrt)
  280. FUNC1(cmath_tan, c_tan)
  281. FUNC1(cmath_tanh, c_tanh)
  282.  
  283.  
  284. static struct methodlist cmath_methods[] = {
  285.     {"acos", cmath_acos, 1},
  286.     {"acosh", cmath_acosh, 1},
  287.     {"asin", cmath_asin, 1},
  288.     {"asinh", cmath_asinh, 1},
  289.     {"atan", cmath_atan, 1},
  290.     {"atanh", cmath_atanh, 1},
  291.     {"cos", cmath_cos, 1},
  292.     {"cosh", cmath_cosh, 1},
  293.     {"exp", cmath_exp, 1},
  294.     {"log", cmath_log, 1},
  295.     {"log10", cmath_log10, 1},
  296.     {"sin", cmath_sin, 1},
  297.     {"sinh", cmath_sinh, 1},
  298.     {"sqrt", cmath_sqrt, 1},
  299.     {"tan", cmath_tan, 1},
  300.     {"tanh", cmath_tanh, 1},
  301.     {NULL,        NULL}        /* sentinel */
  302. };
  303.  
  304. void
  305. initcmath()
  306. {
  307.     object *m, *d, *v;
  308.     
  309.     m = Py_InitModule("cmath", cmath_methods);
  310.     d = getmoduledict(m);
  311.     dictinsert(d, "pi", v = newfloatobject(atan(1.0) * 4.0));
  312.     DECREF(v);
  313.     dictinsert(d, "e", v = newfloatobject(exp(1.0)));
  314.     DECREF(v);
  315. }
  316.